/*************************************************************
*
* M3D_MATH.C
*
* Simple math functions
* Solution of quadratic equations
* Transformation of the angles to different ranges
*************************************************************/

#include <stdafx.h>

//#include <stdlib.h>
//#include <limits.h>
//
//#include <ap_comm.h>
//#include <m3d.h>
//#include <fileio.h>

m3dRND_GEN m3dRndGen;

/*************************************************************
* m3dRND_GEN class members
*************************************************************/
/*************************************************************
* m3dRND_GEN::Next ()
*
* generates next random seed
* IN  : none
* OUT : none
* NOTE: this version is compatible with Microsoft random
*       generator
*************************************************************/
ULONG m3dRND_GEN::Next (void)
{
   ULONG _seed;
   _seed =  seed * 214013UL + 2531011UL;
    seed = _seed;
   return _seed ;
}
/*************************************************************
* m3dRND_GEN::RndInt ()
*
* generates next seed and returns random integer from range [0 .. RAND_MAX]
* IN  : none
* OUT : none
* RET : random integer from range [0 .. RAND_MAX]
*************************************************************/
int m3dRND_GEN::RndInt (void)
{
   return (Next() >> 16) & M3D_RAND_MAX;
}

/******************************************************************************/
/* m3dRND_GEN::RndIntNew
/*
/******************************************************************************/
int m3dRND_GEN::RndIntNew(int prev)
{
   int res;
   do {
      res = RndInt();
   } while (res == prev);
   return res;
}

/*************************************************************
* m3dRND_GEN::RndIntMax ()
*
* generates next seed and returns random integer from range [0 .. max]
* IN  : none
* OUT : none
* RET : random integer from range [0 .. max]
*************************************************************/
int m3dRND_GEN::RndIntMax (int max)
{
   ASSERT(max < M3D_RAND_MAX);
   ULONG mul = (max+1) * ((Next() >> 16) & M3D_RAND_MAX);
   return min(max, int(mul / M3D_RAND_MAX));
}
/*************************************************************
* m3dRND_GEN::RndIntRange ()
*
* generates next seed and returns random float from range [a, b]
* IN  : none
* OUT : none
* RET : random int from range [a, b]
*************************************************************/
int m3dRND_GEN::RndIntRange (int a, int b)
{
   return a + RndIntMax(b-a);
}

/*************************************************************
* m3dRND_GEN::RndFloat ()
*
* generates next seed and returns random float from rang [0, 1]
* IN  : none
* OUT : none
* RET : random float from range [0, 1]
*************************************************************/
float m3dRND_GEN::RndFloat (void)
{
   return (float)((double)Next() / UINT_MAX);
}
/*************************************************************
* m3dRND_GEN::RndFloatMax ()
*
* generates next seed and returns random float from rang [0, max]
* IN  : none
* OUT : none
* RET : random float from rang [0, 1]
*************************************************************/
float m3dRND_GEN::RndFloatMax (float max)
{
   return (float)((double)Next() * max / UINT_MAX);
}
/*************************************************************
* m3dRND_GEN::_RndFloatMax ()
*
* generates next seed and returns random float from rang [-max, max]
* IN  : none
* OUT : none
* RET : random float from rang [0, 1]
*************************************************************/
float m3dRND_GEN::_RndFloatMax (float max)
{
   return -max + 2.f * (float)((double)Next() * max / UINT_MAX);
}

/*************************************************************
* m3dRND_GEN::RndFloatRange ()
*
* generates next seed and returns random float from rang [0, max]
* IN  : none
* OUT : none
* RET : random float from rang [0, 1]
*************************************************************/
float m3dRND_GEN::RndFloatRange (float a, float b)
{
   return a + (float)((double)Next() * (b-a) / UINT_MAX);
}


/*************************************************************
* m3dRND_GEN::IsProbability()
*
* Get a random value in [0..100] and return TRUE if it is
* less than threshold
*************************************************************/
BOOL m3dRND_GEN::IsProbability (float threshold)
{
   if (RndFloat()*100.f < threshold) {
      return TRUE;
   }
   return FALSE;
}

/******************************************************************************
* m3dRND_GEN::RndBasicCoord ()
*
* gets coefficients for a random point in a triangle
* OUT: pA - pointer to first coefficient
*      pB - pointer to second coefficient
******************************************************************************/
void m3dRND_GEN::RndBasicCoord (float *pA, float *pB)
{
   float a, b;

   ASSERT(pA != NULL);
   ASSERT(pB != NULL);

   a = RndFloat();
   b = RndFloat();
   if (a + b <= 1.f) {
      *pA = a;
      *pB = b;
   } else {
      *pA = 1.f - b;
      *pB = 1.f - a;
   }
   return;
}
/*************************************************************
* m3dRND_GEN::RndVUnit ()
*
* OUT: v - a random unit vector
*
* Make a random unit vector
*************************************************************/
void m3dRND_GEN::RndVUnit (m3dV *v)
{
   v->x =  RndFloat() - 0.5f; 
   v->y =  RndFloat() - 0.5f;        
   v->z =  RndFloat() - 0.5f; 
   m3dNormalize(v);
   return;
}
/*************************************************************
* m3dRND_GEN::RndVUnitConus ()
* 
* IN:   vBase - source vector
*       angle - restrict angle
* OUT:  vRand - random vector that belongs to conus
*
* Generate random unit vector that belongs to conus
*************************************************************/
void m3dRND_GEN::RndVUnitConus (const m3dV *vBase, float angle, m3dV *vRand)
{
   m3dV   vBaseScale, vTangScale;
   m3dV   v, vTang;
   float  a;

   if (m3dLengthVector(vBase) < M3D_EPSILON) {
      m3dMakeVZero(vRand);
      return;
   }
   RndVUnit(&v);
   m3dTangentComponent(&v, vBase, &vTang);
   if (m3dLengthVector(&vTang) > M3D_EPSILON)   {
      m3dNormalize(&vTang);
      a = m3dDeg2Rad(RndFloatMax(angle));
      m3dScaleVector(vBase,  m3dCos(a), &vBaseScale);
      m3dScaleVector(&vTang, m3dSin(a), &vTangScale);
      m3dAddVector(&vBaseScale, &vTangScale, vRand);
      m3dNormalize(vRand);
   } else   {
      *vRand = *vBase;
   }

   return;
}
/*************************************************************
* m3dRND_GEN::RndVTang ()
*    Generate random vector perprndicular to normal
* IN:  norm  - normal vector
* OUT: vRand - random vector that belong conus
* RET: none
*************************************************************/
void m3dRND_GEN::RndVTang (m3dV *norm, m3dV *randDir)
{
   m3dV normDir = *norm;
   
   if (!_m3dNormalize(&normDir)) {
      normDir = m3dVUnitY;
   }
   // generate direction perpendicular to dir vector
   RndVUnit(randDir);
   m3dTangentComponent(randDir, &normDir, randDir);
   m3dNormalize(randDir);
   return;
}

/*************************************************************
* m3dRND_GEN::RndVCylinder()
*    Generate random vector within a given cylinder
* IN:  vBase   - cylinder axis
*      radius  - cylinder radius
*      height  - cylinder height
* OUT: vRand   - random vector that belong it
* RET: none
*************************************************************/
void m3dRND_GEN::RndVCylinder(const m3dV *vpBase, float radius, float height, m3dV *vpRand)
{
   m3dV  vTemp;

   ASSERT(vpBase);
   ASSERT(radius > 0);
   ASSERT(height > 0);

   // check if base vector is zero
   if (m3dLengthVector(vpBase) < M3D_EPSILON) {
      m3dMakeVZero(vpRand);
      return;
   }

   // generate direction perpendicular to dir vector
   RndVUnit(&vTemp);
   m3dTangentComponent(&vTemp, vpBase, &vTemp);
   m3dNormalize(&vTemp);
   _m3dScaleVector(&vTemp, RndFloatMax(radius));

   // generate random height vector
   m3dScaleVector(vpBase, RndFloatMax(height), vpRand);

   // prepare result vector
   m3dAddVector(&vTemp, vpRand, vpRand);

   return ;
}

/*************************************************************
* m3dRND_GEN::RndTriPoint()
*    Generate random point in triangle
* IN:  tri   - triangle
* OUT: vRand - random vector that belong it
* RET: none
*************************************************************/
void m3dRND_GEN::RndTriPoint(m2dTRI *tri, m2dV *vRand)
{
   m2dV  vA, vB;
   float kA, kB;
   
   // generate random position for the particle
   RndBasicCoord(&kA, &kB);
   m2dSubtractVector(&tri->v2, &tri->v1, &vA);
   m2dSubtractVector(&tri->v3, &tri->v1, &vB);
   m2dLineCombine(&vA, kA, &vB, kB, vRand);
   _m2dAddVector(vRand, &tri->v1);
   return;
}

/*************************************************************
* m3dRND_GEN::RndTriPoint()
*    Generate random point in triangle
* IN:  tri   - triangle
* OUT: vRand - random vector that belong it
* RET: none
*************************************************************/
void m3dRND_GEN::RndTriPoint(m3dTRI *tri, m3dV *vRand)
{
   m3dV  vA, vB;
   float kA, kB;
   
   // generate random position for the particle
   RndBasicCoord(&kA, &kB);
   m3dSubtractVector(&tri->v2, &tri->v1, &vA);
   m3dSubtractVector(&tri->v3, &tri->v1, &vB);
   m3dLineCombine(&vA, kA, &vB, kB, vRand);
   _m3dAddVector(vRand, &tri->v1);
   return;
}



/*************************************************************
* m3dRND_DISTR class members
*************************************************************/
/*************************************************************
* m3dRND_DISTR::SetWeigth ()
*    Set weight of the item in the distribution
*************************************************************/
//void m3dRND_DISTR::SetWeigth (int weigthNmb, float weigth)
//{
//   SetKpArg(weigthNmb + 1, weigth);
//   return;
//}
///*************************************************************
//* m3dRND_DISTR::UnifyWeigths ()
//*    Define number of weights and unify them
//*************************************************************/
//void m3dRND_DISTR::UnifyWeigths(void)
//{
//   /*
//   if (_nWeights == -1) {
//      nWeights = NKp() - 1;
//   } else {
//      nWeights = _nWeights;
//   }
//   SetNKp(nWeights + 1);
//   */
//   int nWeights = NKp() - 1;
//   
//   // prepare a spline wich has keypoint argument distribution
//   // according to the distribution of item weights
//   // in [0...1] interval
//   SetKpArg(0, 0.0f);
//   for (int i = 0; i < nWeights; i ++) {
//      SetKpArg(i + 1, GetKpArg(i) + GetKpArg(i + 1));
//   }
//   for (int i = 0; i < nWeights; i ++) {
//      SetKpArg(i, GetKpArg(i) / GetKpArg(nWeights));
//   }
//   SetKpArg(nWeights, 1.0f);
//
//   return;
//}
///*************************************************************
//* m3dRND_DISTR::GetSample ()
//*    Get random sample
//*************************************************************/
//void m3dRND_DISTR::GetSample (void *sample)
//{
//   int  ind;
//
//   ind = FindSeg(m3dRandom());
//   switch(ind) {
//      case FIND_SEG_LEFT:
//         ind = 0;
//         break;
//      case FIND_SEG_RIGHT:
//         ind = NKp() - 2;
//         break;
//   }
//   *((int*)sample) = ind;
//   return;
//}
//
//
///*************************************************************
//* m3dRND_DISTR_SURF class members
//*************************************************************/
///*************************************************************
//* m3dRND_DISTR_SURF::SetWeigth ()
//*    Set weight of the item in the distribution
//*************************************************************/
//void m3dRND_DISTR_SURF::SetWeigth (int weigthNmb, m3dTRI *pTri)
//{
//   SetKpValue(weigthNmb, pTri);
//   return;
//}
/*************************************************************
* m3dRND_DISTR_SURF::SetWeigth ()
*    Set weight of the item in the distribution
*************************************************************/
//void m3dRND_DISTR_SURF::SetWeigth (int weigthNmb, m3dV *v1, m3dV *v2, m3dV *v3)
//{
//   m3dTRI   tri;
//   tri.v1 = *v1;
//   tri.v2 = *v2;
//   tri.v3 = *v3;
//   SetKpValue(weigthNmb, &tri);
//   return;
//}
///*************************************************************
//* m3dRND_DISTR_SURF::UnifyWeigths ()
//*    Define number of weights and unify them
//*************************************************************/
//void m3dRND_DISTR_SURF::UnifyWeigths(void)
//{
//   /*
//   if (_nWeights == -1) {
//      _nWeights = NKp() - 1;
//   } 
//   */
//   for (int i = 0; i < NKp() - 1; i ++) {
//      m3dTRI tri;
//      GetKpValue(i, &tri);
//      m3dRND_DISTR::SetWeigth(i, m3dCalcTriArea(&tri.v1, &tri.v2, &tri.v3));
//   }
//   m3dRND_DISTR::UnifyWeigths();
//   return;
//}
///*************************************************************
//* m3dRND_DISTR_SURF::GetSample ()
//*    Get random sample
//*************************************************************/
//void m3dRND_DISTR_SURF::GetSample (void *sample)
//{
//   m3dRND_DISTR_SURF_SAMPLE *pSample = (m3dRND_DISTR_SURF_SAMPLE*)sample;
//
//   m3dRND_DISTR::GetSample(&pSample->ind);
//
//   m3dMakeTriPointRandom((m3dTRI*)&GetRaw(pSample->ind).value, &pSample->point);
//   
//   return;
//}

/*************************************************************
'* m3dRND_EXCLUSIVE::Init()
'* 
'*************************************************************/
void m3dRND_EXCLUSIVE::Init(int size)
{
   // decide by size
   if (size > 1) {// have some
      topIdx  = size - 1;
      lastIdx = m3dRandIntMax(topIdx);
   } else { // too few => no selection
      lastIdx = topIdx  = 0;
   }
}

/*************************************************************
'* m3dRND_EXCLUSIVE::Get()
'* 
'* RET : new index within [0, size-1]
'*************************************************************/
int  m3dRND_EXCLUSIVE::Get(void)
{
   // determined or not
   switch (topIdx) {
      case 0: break;// no choice no worry

      case 1: lastIdx = 1-lastIdx; break; // flip-flop

      default: // now we have a selection
         int shift = 1 + m3dRandIntMax (topIdx - 1);
         lastIdx = (lastIdx + shift) % (topIdx + 1);
   }

   return lastIdx;
}

DWORD m3dVRand(const m3dV &vv, DWORD stride) 
{
   float intPart;
   m3dV v = vv*10;
   v.Set((float)m3dIntFloor(v.x), (float)m3dIntFloor(v.y), (float)m3dIntFloor(v.z));
   v /= 1000;
   DWORD r1 = (DWORD)(modf(v.x, &intPart) * 65535.f * 9241);
   DWORD r2 = (DWORD)(modf(v.y, &intPart) * 65535.f * 7057);

   DWORD r3 = (DWORD)(modf(v.z, &intPart) * 65535.f * 11719);
   DWORD r = r1 ^ (r2 << 8) ^ (r3 << 16) ^ (r2 >> 7) ^ (r3 >> 3);
   return r * stride;
}
int m3dVRandIntRange(const m3dV &v, int a, int b, DWORD stride)
{
   int d = b - a;
   DWORD r = m3dVRand(v, stride) & M3D_RAND_MAX;
   DWORD mul = d * r;
   return a + mul / M3D_RAND_MAX;
}

float m3dVRandRange(const m3dV &v, float a, float b, DWORD stride)
{
   DWORD r = m3dVRand(v, stride);
   return a + (float)((double)r * (b-a) / UINT_MAX);
}
void  m3dMakeVUnitVRandom(const m3dV &cc, m3dV *v, DWORD stride) 
{
   m3dV c = cc;
   v->x = m3dVRandRange(c, -0.5f, 0.5f, stride); c.Set(c.z, c.x, c.y);
   v->y = m3dVRandRange(c, -0.5f, 0.5f, stride); c.Set(c.z, c.x, c.y);
   v->z = m3dVRandRange(c, -0.5f, 0.5f, stride); 
   m3dNormalize(v);
}  

void m3dMakeVUnitConusVRandom(const m3dV &cc, m3dV *vBase, float angle, m3dV *vRand, DWORD stride)
{
   m3dV   vBaseScale, vTangScale;
   m3dV   v, vTang;
   float  a;
   if (m3dLengthVector(vBase) < M3D_EPSILON) {
      m3dMakeVZero(vRand);
      return;
   }
   m3dV c = cc;
   m3dMakeVUnitVRandom(c, &v, stride);
   m3dTangentComponent(&v, vBase, &vTang);
   if (m3dLengthVector(&vTang) > M3D_EPSILON)   {
      m3dNormalize(&vTang);
      c.Set(c.y, c.x, c.z);
      a = m3dDeg2Rad(m3dVRandRange(c, 0, angle, stride));
      m3dScaleVector(vBase,  m3dCos(a), &vBaseScale);
      m3dScaleVector(&vTang, m3dSin(a), &vTangScale);
      m3dAddVector(&vBaseScale, &vTangScale, vRand);
      m3dNormalize(vRand);
   } else   {
      *vRand = *vBase;
   }
}
//
// End-of-file M3D_RAND.CPP
//
